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We present the exact Bethe ansatz solution for the two-dimensional BCS pairing Hamiltonian 
with p x + ip y symmetry. Using both mean- field theory and the exact solution we obtain the ground- 
state phase diagram parameterized by the filling fraction and the coupling constant. It consists of 
three phases denoted weak coupling BCS, weak pairing, and strong pairing. The first two phases are 
separated by a topologically protected line where the exact ground state is given by the Moore- Read 
pfaffian state. In the thermodynamic limit the ground-state energy is discontinuous on this line. The 
other two phases are separated by the critical line, also topologically protected, previously found 
by Read and Green. We establish a duality relation between the weak and strong pairing phases, 
whereby ground states of the weak phase are "dressed" versions of the ground states of the strong 
phase by zero energy (Moore-Read) pairs and characterized by a topological order parameter. 
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In 1957, Bardeen, Cooper and Schrieffer published an 
epoch defining paper giving a microscopic explanation of 
the properties of superconducting metals at low tempera- 
tures [l( . The model was based on a reduced Hamiltonian 
which describes the pairing interaction between conduc- 
tion electrons. The original study of the BCS model was 
formulated in the grand canonical ensemble and solved 
with a mean-field approximation. In 1963 Richardson 
derived the exact solution of the reduced BCS Hamil- 
tonian with s-wave symmetry in the canonical ensemble 
[2] ■ This solution was largely unnoticed until its rediscov- 
ery in the theoretical studies of ultrasmall metallic grains 
in the 1990's, where it was employed to understand the 
crossover between the fluctuation dominated regime and 
the fully developed superconducting regime (for a review 
sec [3]). The exact solution for the s-wave BCS model 
is related to the Gaudin spin Hamiltonians and their in- 
tegrability can be understood in the general framework 
of the Quantum Inverse Scattering Method 0, Q . These 
later developments allowed for an exact computation of 
various correlators @, H, @] , and led to generalizations of 
the Richardson-Gaudin models with applications to con- 
densed matter and nuclear physics [3j, |8| . 

In this Letter we analyze the two-dimensional BCS 
model where the symmetry of the pairing interaction is 
Px+iPy (hereafter referred to as p+ip). The Hamiltonian 
of the model is 
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where Ck, c k are destruction and creation operators of 
2D spinless or polarised fermions with momentum k, m 
is their mass and G is a dimensionless coupling constant 
which is positive for an attractive interaction. The p + ip 
model has attracted considerable attention due to the 



connection with the Moore- Read pfaffian state arising in 
the quantum Hall effect at filling fraction 5/2 [9j, which 
has been proposed to support non-abelian anyons allow- 
ing for topological quantum computation [Ig|, [H| ■ Mo- 
tivated by these considerations, concrete proposals for 
engineering the p + ip form of the pairing interaction 
have been formulated in the context of cold fermi gases 
I2I 13j |. Here we will study the model through the ex- 



act Bethe ansatz solution, which is presented for the first 
time. We remark that exact solvability holds indepen- 
dent of the choice for the ultraviolet cut-off, which we de- 
note as oj, and independent of the distribution of the mo- 
menta k. In particular this means that a one-dimensional 
system is obtained by simply setting all k y — 0. Unless 
stated otherwise, all discussions below deal with finite 
particle numbers in a finite-sized system. 

Using the standard mean-field theory approach Read 
and Green have shown the existence of a second-order 
jhase transition governed by the chemical potential /i 
Til ]. Adopting the terminology of [3], this transition 
takes place between a weak pairing phase (/j, > 0), the 
ground state (GS) of which behaves as the Moore-Read 
pfaffian state at long distances, and a strong pairing 
phase (for /1 < 0). The spectrum of Bogolioubov quasi- 
particles is gapless at fj, = 0. The GS of the weak pair- 
ing phase also has a non-trivial topological structure in 
k-space, as shown by Volovik [l5| . However in the mean- 
field analysis the weak pairing GS is continuously con- 
nected to the weak coupling BCS GS [14j . 

Our goal is to re-examine the properties of the p + ip 
model. Through this study we will achieve the following: 
i) From the mean-field results the ground-state phase di- 
agram will be determined, comprising of the weak cou- 
pling BCS, weak pairing, and strong pairing phases; ii) 
A duality relation between the weak pairing and strong 
pairing phases will be shown to exist; iii) From the Bethe 
ansatz solution the duality will be formulated in terms 
of a dressing relation involving zero energy Moore-Read 
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(MR) pairs; iv) Dressing of the vacuum will be seen to 
give the boundary line between weak coupling BCS and 
weak pairing phases, representing a zeroth-order quan- 
tum phase transition when the thermodynamic limit is 
taken (cf. [l6| for analogous zeroth-order thermal phase 
transitions); iv) The weak pairing phase will be shown 
to have a non-trivial topological structure, related to the 
dressing operation, which will be quantified by a winding 
number. 

Before presenting the exact solution of the Hamiltonian 
we first extend the mean- field results reported in 
The BCS order parameter associated to ([T]) is 
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in terms of which the Hamiltonian ([I]) can be approxi- 
mated as (up to an additive constant) 
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where £k = k 2 /2m — /i/2 and /x/2 is the chemical po- 
tential. This Hamiltonian can be diagonalized by a Bo- 
goliubov transformation. The gap A = |A| and chemical 
potential are the solutions of the equations 
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where we have set m — 1, Lis the total number of energy 
levels and M is the number of Cooper pairs. The set K + 
denotes the set of momenta where k x > and any k v , so 
that we avoid overcounting of energy levels. The mean- 
field expression for the GS energy is (accounting for the 
constant term missing in ([3])) 
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Projection of the grand-canonical GS wave function onto 
a fixed number of M pairs gives 
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where g{k) = {2E(k) - k 2 + fx)/((k x + 
E(k) is the quasiparticle energy spectrum 
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and 
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Note that the spectrum is gapless at [i = as |k| — > 0. 
Furthermore, the behaviour of g(k) as k — > depends on 
the sign of fi [lj] , 
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FIG. 1: Ground-state phase diagram of the p + ip model 
in terms of the inverse coupling 1/g and filling fraction x = 
M/L. Phase boundaries are given by the Read-Green line 
(fj. = 0) and the Moore-Read line (/x = A 2 /4) ^fj]. The phase 
boundaries are independent of the choice of the momentum 
distribution, and independent of the ultraviolet cut-off. 



In real space the state (J7J takes the form of a pfaffian 
V>(ri, . . . ,r 2 Af) = A\g(ri - r 2 ) . ..g(r 2 M-i -tzm)] (10) 

where A denotes the antisymmetrization of the positions 
and g(r) is the Fourier transform of <?(k). We will refer 
to the case /i = as the Read-Green (RG) state. For 
fi > the large distance behaviour is g(r) ~ l/(x + iy), 
which asymptotically reproduces the MR state [14| . 

The solution of the equations (|4I5[) for L energy levels 
and M number of pairs can be classified, with the corre- 
sponding phase diagram given in fig. [T] parameterized by 
the filling fraction x = M/L and the rescaled coupling 
constant g — GL. We now demonstrate how the topo- 
logical aspects of the phase diagram can be deduced in a 
transparent manner. From ([5|) we see that /i — imposes 
the relation xrg = (1 — 5 _1 )/2- This result is completely 
independent of the momentum distribution and choice of 
cut-off, reflecting the topological nature of the transition 
discussed in i.e., the boundary line is protected from 
perturbations of the system which alter the distribution 
of the momenta. Furthermore we identify a second topo- 
logical boundary by setting fi = A 2 /4, which from (JSJ) 
gives E = 0, again independent of the momenta. Using 
(|4l5p it is found this occurs when xmr = 1 — g 1 - Later 
we will show that in this instance the GS is a discrete 
analogue of the MR state mentioned earlier, which in the 
thermodynamic limit is exactly the MR state. 

A further notable k-independent property of the phase 
diagram is the existence of a "duality" between a point 
(g,xi) in the weak pairing regime and another point 
(g,xu) in the strong pairing regime related by 
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which necessarily can only hold for rational values of g. 
In the mean-field analysis this duality means that the 
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corresponding solutions are related by fii = —fJ-ii and 
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■ 2 fin, such that the GS energies satisfy 



Ei = En according to ([6]). The RG state is self-dual, 
whereas the MR state is dual to the vacuum. This duality 
is apparent in the exact solution where it will be shown 
to be related to a dressing operation mentioned in the 
introduction. 

The detailed derivation of the exact Bethe ansatz so- 
lution will be presented elsewhere. Here we simply men- 
tion that the technical aspects follow the derivation of 
the s-wave model solution through the Quantum Inverse 
Scattering Method, as described in @, Q . The only fun- 
damental difference is that the i?-matrix solution of the 
Yang-Baxter equation used to solve the p + ip model is 
the trigonometric XXZ solution, in contrast to the ra- 
tional XXX solution used for the s-wave model. 

We again set m = 1. The exact eigenstates of the 
Hamiltonian with M fermion pairs are given by 



M 



J\C(y 3 )\0), C(y) 



E 



yj J 



k 2 



y 



(12) 



where the rapidities yj, j — 1,...,M satisfy the Bethe 
ansatz equations (BAE) 
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with 2q = l/G — L + 2M — 1. The total energy of the 
state (fT2"j) is given by 
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Numerical solutions of the Bethe ansatz equations indi- 
cate that there are no unpaired fermions in the GS when 
the fermion number is even. In fig. [3 we present numer- 
ical GS solution of eqs. (fTB")) . This solution is obtained 
starting from the initial condition yj — > (1 + G)k 2 (j = 
1, . . . , M) as G — > 0, with the k chosen to fill the Fermi 
sea. As g increases, the roots yj closest to the Fermi level 
become complex pairs. When g approaches the MR-line 
the roots bend towards the origin (as shown in fig. [2j) 
and at the value g~ x = 1 — x all the roots collapse onto 
the origin (not shown). At larger values of g one enters 
the weak pairing phase where all the roots are non-zero, 
except at some rational values of g where a fraction of 
the roots collapse again. Finally, in the strong pairing 
regime all the roots become real and they belong to an 
interval on the negative real axis. 

Looking closer at the weak pairing phase, one can check 
that the M\y roots yj can be split into Mq vanishing roots 
and Ms non-zero roots provided 
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FIG. 2: Numerical solutions for the ground-state roots 
yj (j = 1, . . . , M) of the BAEs {TJJ with M = 31, L = 62, 
and < g < 1.99. For this range of couplings the system is in 
the weak coupling BCS phase. The maximal pairing energy 
(the cut-off) is denoted u. At the critical coupling g = 2 (MR 
line) all the roots collapse to the origin (not shown). 



Moreover the Ms non-zero roots satisfy the BAE (fTB"]) in 
the strong pairing region. Altogether this implies that 
given an eigenstate, say \S), in the strong pairing regime 
then one can dress it with Mq MR pairs (as given by 
(fTS")) ) obtaining an eigenstate \W) in the weak pairing 
phase with the same energy, i.e. 

H\S) = E\S) H\W) = H[C{{))] Ma \S) = E\W). 

Noticing that the filling fraction of the strong pairing 
state \s xs = Ms/L and that of the weak pairing state is 
xw = (Mq+Ms)/L, we find that eq. (fH)]) coincides with 
the duality relation (jlip . The physical picture we obtain 
from this discussion is that the fermion pairs forming the 
GS in the weak pairing phase are of two types: strong 
localized pairs with negative energy and the delocalized 
MR pairs with zero energy. This picture is substantially 
different from the projected mean- field wavefunction ((7|, 
which is more akin to a condensate of Cooper pairs in the 
same one-particle state. An exception to this occurs on 
the MR line, where the projected mean-field and exact 
wavefunctions are identical. We see from (TT2"]) that when 
all roots of the Bethe ansatz equations are zero, the GS 
is a discrete analogue of the MR state with zero energy 
in agreement with mean-field theory. 

We reiterate that until now all our analysis has been 
in the context of finite-sized systems, and in particular 
the topological (i.e. k-independent) nature of the duality 
(jlljl is not dependent on taking the thermodynamic limit. 
In going to the thermodynamic limit we take L, M — > oo, 
G — * with x = Mj L and g = GL fixed. A peculiar fea- 
ture of the MR line is the discontinuity of the GS energy 
E(g, x) in the thermodynamic limit as the filling frac- 
tion x approaches the value Xmr from the weak pairing 
region. To derive this result, for finite L we take the one- 
pair state and dress it to give the dual GS in the weak 
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pairing region. The filling xj of the dressed state is given 
by fTTjh setting xjj — l/L, i.e. xi = xmr —1/L, which 
approaches x^r = 1 — l/s as L — > 00. Since the MR 
pairs carry no energy, the GS energy of the dressed state 
coincides with the one-pair energy. To compute this en- 
ergy we consider the BAE for one Cooper pair and take 
the continuum limit (i.e. eq. (fl3|) with M = 1). Settting 
e = k 2 and u> as the cut-off, for simplicity we take the 
momentum distribution to be that for free particles in 
two dimensions, i.e. p(e) = u)~ x . This leads to 
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This equation has a unique negative energy solution y < 
satisfying 



1 - - = - log 

g w \y-w 



which we denote as y — 8(g). From here one derives 
the aforementioned discontinuity on the MR line Xmr — 

1-9- 1 , 

lim E(g,xi) = £(g) ^ E(g, xmr) = 0, 

L — >oo 

which may be described as a zeroth-order quantum phase 
transition. To the best of our knowledge, this is the first 
example of a zeroth-order quantum phase transition in a 
many-body system. We have also numerically analyzed 
the excited states on the MR line obtained by blocking 
the energy levels which are occupied by unpaired elec- 
trons. These excitations have a gap whose value agrees 
with the mean-field result, suggesting that only the RG 
line is gapless, consistent with mean-field theory predic- 
tions. 

As mentioned in the introduction, the mean-field solu- 
tion shows that the weak pairing phase has a non-trivial 
topological structure in k-space LJ, [15J . This structure 
can be characterized by the winding number w of the 
mean- field wavefunction <?(k) = g x (k) + ig y ()s.), and it is 
given by, 
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One finds that w — for /1 < (i.e. strong pairing 
phase), while w = +1 for /z > (i.e. weak pairing and 
weak coupling BCS phases) [3, EH- The existence of an 
exact solution of the model calls for a generalization of w 
applicable to the many-body wavefunction of the model 
■0(ki, . . . , 1cm), where kj (i = 1, . . . , M) are the distinct 
momenta of the pairs. This generalization consists in 
replacing g(k) in IT6|) by ^(k + Ci, . . . , k + cm), where 
Cj ^ ci Vj, I are a set of distinct constants. With this def- 
inition we find that w vanishes for the exact ground-state 
wavefunction except in the weak pairing region where it 



coincides with the number of MR pairs. Hence w pro- 
vides a non-trivial topological order parameter for the 
weak pairing phase which is zero in the other two phases. 

In summary, we have provided the exact Bethe ansatz 
solution for the BCS model with p + ip pairing. Us- 
ing this we have investigated the ground-state phase dia- 
gram, whose structure is richer than previously supposed. 
We have found that the weak pairing region is dual to 
the strong pairing region, the duality being encoded in 
a dressing transformation between GS of the two phases 
by means of zero energy MR pairs. The MR state ob- 
tained by dressing the vacuum is the exact GS on a line 
in the phase diagram. The MR line separates the weak 
pairing and weak coupling BCS regions, and while the 
gap does not vanish on it, the GS energy is discontinu- 
ous in the thermodynamic limit. We have also found a 
topological order parameter that characterizes the weak 
pairing phase. An important future issue is to explore 
how vortices (e.g. see [lOj) can be incorporated into a 
similar model to the one studied in this Letter. 
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